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Abstract 

A statistically principled way of conducting weighted networJj analysis is still lacking. 
Comparison of different populations of weighted networks is hard because topology is inher- 
ently dependent on wiring cost, where cost is defined as the number of edges in an unweighted 
graph. In this paper, we evaluate the benefits and limitations associated with using cost- 
integrated topological metrics. Our focus is on comparing populations of weighted undirected 
graphs using global efficiency. We evaluate different approaches to the comparison of weighted 
networks that differ in mean association weight. Our key result shows that integrating over 
cost is equivalent to controlling for any monotonic transformation of the weight set of a 
weighted graph. That is, when integrating over cost, we eliminate the differences in topology 
that may be due to a monotonic transformation of the weight set. Our result holds for any 
unweighted topological measure. Cost-integration is therefore helpful in disentangling differ- 
ences in cost from differences in topology. By contrast, we show that the use of the weighted 
version of a topological metric does not constitute a valid approach to this problem. Indeed, 
we prove that, under mild conditions, the use of the weighted version of global efficiency is 
equivalent to simply comparing weighted costs. Thus, we recommend the reporting of (i) 
differences in weighted costs and (ii) differences in cost-integrated topological measures. We 
demonstrate the application of these techniques in a re-analysis of an fMRI working memory 
task. Wc also provide a Monte Carlo method for approximating cost-integrated topological 
measures. Finally, we discuss the limitations of integrating topology over cost, which may pose 
problems when some weights are zero, when multiplicities exist in the ranks of the weights, 
and when one expects subtle cost-dependent topological differences, which could be masked 
by cost-integration. 

KEYWORDS: Connectivity, Correlation matrix, Cost-integration, Global Efficiency, Monte Carlo 
integration, Networks, Small-world. 

1 Introduction 

In the last decade, the biological and physical sciences have witnessed a proliferation of publications 
adopting a network approach to a wide range of questions. This interest in networks was originally 
stimulated by the seminal works of Watts and Strogatz (1998) and Barabasi and Albert (1999), 
who introduced the concepts of small-world and scale-free networks, respectively. Some of these 
ideas have been adopted in neuroscience at both a theoretical (Sporns et al., 2000, 2004) and 
experimental level (Eguiluz et al., 2005). Most of the research in this area has attempted to 
classify the topology of brain networks based on anatomical or functional data (see, for example, 
Achard et al., 2006, Achard and Bullmore, 2007, He et al., 2007). 

A question that naturally arises from such applications of graph theory is whether or not the 
topological properties of these brain networks are stable across different populations of subjects 
or across different cognitive and behavioral tasks. A common hypothesis that neuroscientists 
may wish to test is whether the small-world properties of a given brain network are conserved 
when comparing patients and healthy controls. Bassett et al. (2008), for example, have studied 
differences in anatomical brain networks between healthy controls and patients with schizophrenia. 
Other authors have evaluated whether the topological properties of functional networks vary with 
different behavioral tasks (van den Heuvel et al., 2009, De Vico Fallani et al., 2008, Cecchi et al., 
2007, Astolfi et al., 2009). The properties of brain network topology have also been studied at 
different spatial scales (Bassett et al., 2006) and using different modalities, such as EEC (Pachou 
et al., 2008, Salvador et al., 2008), and fMRI (Achard et al., 2006, Achard and Bullmore, 2007). 
There is therefore considerable interest in comparing populations of networks -which may represent 
different groups of subjects, several conditions of an experiment, or the use of different levels of 
spatial or temporal resolution. We note that such research questions are more likely to arise when 
subject-specific networks can be directly constructed. This has been done in the context of both 
functional and structural MRI (Hagmann et al., 2008, Gong et al., 2009). 

The possibility of conducting rigorous statistical comparison of several populations of networks, 
however, has been hindered by a series of methodological issues, which have not been hitherto sat- 
isfactorily resolved. When considering the question of comparing several populations of networks, 
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Figure 1. How can one disentangle differences in connectivity strength from differences in topology? In 
panel (a), two correlation matrices for two weighted networks differ in their average correlation strengths. 
In panel (b), the same correlation matrices have been thresholded at the same value, producing graphs 
with different cost levels. In all matrices, black indicates null values, and white denotes entries equal to 
unity. 



two main problems arise. Firstly, we are faced with the inherent intertwining of connectivity 
strength (i.e. wiring cost) with network topology. Most topological metrics used to compare net- 
works are sensitive to differences in these graphs' number of edges. Drawing comparisons on 
the sole basis of topology therefore requires some level of control of cost discrepancies between 
these network populations. Secondly, this issue is compounded by the fundamental division be- 
tween weighted and unweighted graphs. The problem of disentangling differences in connectivity 
strength from topological differences therefore needs to be resolved in a distinct manner depending 
on whether weighted or unweighted graphs are being considered. The focus, in this paper, will be 
on weighted networks since these are more likely to be found in the biomedical sciences than their 
unweighted counterparts. 

Historically, however, network analyses have concentrated on unweighted graphs. The applica- 
tion of graph theory to biological and artificial ncitworks was originally motivated by the discrete 
nature of the problems of interest. Both Watts and Strogatz (1998) and Barabasi and Albert 
(1999) mainly considered binary relations between sets of elements, which readily produced ad- 
jacency matrices that could then be used to construct unweighted graphs. Watts and Strogatz 
(1998) matched some networks of interest with their random and regular equivalents. In their case, 
the matching procedure ensured that both random and regular networks possessed the same total 
number of nodes and edges as the original graph. Current practice in MRI-based neuroscience and 
other biomedical applications, however, tends to produce weighted connectivity networks. This is 
because MRI data take values on a continuous scale, which lends itself to the application of real- 
valued measures of association, such as the correlation coefficient or the synchronization likelihood 
among others. While different populations of unweighted networks can readily be compared by 
matching each network with a random network possessing an identical number of edges; there is, 
as yet, no consensus on how to compare populations of weighted networks in a systematic manner. 

This problem can be illustrated with a straightforward example. In panel (a) of Figure 1, 
a pair of weighted networks are represented by their correlation matrices. We are interested in 
comparing the topology of the corresponding weighted graphs. Since these networks differ in 
their mean correlation coefficients, a simple thresholding of these matrices will produce graphs of 
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different levels of cost. Naturally, this is only one of the different thresholding approaches that 
could be adopted. This non-uniqueness arises since graph topology is expressed in the language 
of discrete mathematics, whereas correlation coefficients are real-valued functions. That is, one 
cannot directly adopt concepts originally developed for unweighted graphs for the analysis of 
weighted graphs. 

In this paper, we consider two main approaches to the problem of weighted network comparison. 
Firstly, following other authors, we evaluate the use of weighted topological metrics, which are 
weighted equivalents of graph-theoretical metrics for unweighted networks (Latora and Marchiori, 
2001, Rubinov and Sporns, 2010). Secondly, we consider the utilization of cost-integrated measures 
of topology, where all the possible wiring costs of a network are taken into account. When 
unweighted, a graph's wiring cost is defined as its number of edges. Integrating over wiring 
cost can here be interpreted in statistical terms, as an analog to the Bayesian integration of 
nuisance^ parameters. Doing so, we are averaging out the 'uncertainty' in the choice of a particular 
level of cost. This second family of measures has probably been the most popular to date in 
the neuroscientific literature (see Achard et al., 2006, He et al., 2009b, for instance). However, 
different authors have chosen different integration intervals. We therefore explore the consequences 
of integrating a topological metric with respect to different subsets of the cost interval. 

Note that our approach substantially differs from the one adopted by van Wijk et al. (2010), 
who proposed several formulas relating cost levels and topological measures, such as the charac- 
teristic path length or the clustering coefficient. Instead, in this paper, we are concerned with 
formally deriving what is the effect of integrating a particular topological measures over cost 
levels, in order to assess whether this is a successful manner of disentangling differences in cost 
from differences in topology. In particular, although van Wijk et al. (2010) reviews several ways 
of controlling for differences in cost, they do not consider cost-integration, per se. This paper 
can therefore be seen as a contribution to the literature on weighted network analysis, where we 
formally clarify the utilization of cost-integration when comparing topological metrics. 

The concept of topology in the context of this paper will be defined in a quantitative manner. 
This should be contrasted with the qualitative definition adopted by previous authors, van Wijk 
et al. (2010), for instance, assume that networks that represent different realizations of the same 
'generative model' should be regarded as topologically identical. Several realizations from an 
Erdos-Rcnyi model with fixed edge probability, for example, share a common generative model and 
therefore can be said to have an identical topology. In practice, however, such a generative model 
is unknown. Thus, we will refer to this type of classification as a topological taxon. A taxonony 
of commonly encountered networks may include the random topology of the Erdos-Rcnyi model, 
the regular lattice and the small-world topology among others. Such a nomenclature is qualitative 
because it relies on discrete categories. By contrast, we wish to adopt a quantitative perspective 
on this problem, whereby topology is operationalized in terms of specific topological properties 
such as the clustering coefficient (CC), for instance. In this perspective, two Erdos-Renyi models 
with identical edge probability may display different levels of global and local efficiencies and 
will therefore be considered to have distinct topological properties. Therefore, we distinguish 
between a qualitative approach based on topological taxonony and a quantitative approach based 
on topological properties. Given that generative models are latent, our quantitative definition of 
topology appears better suited to the empirical study and comparison of complex networks. 

The above definition of network topology, however, assumes that the networks under com- 
parison have identical numbers of vertices and edges. When this is not the case, or when one is 
comparing two popiilations of weighted networks, the question of whether or not these networks 
have similar topological properties become arduous. Our main aim, in this paper, is therefore 
to identify the situations within which one can safely conclude that different weighted networks 
share the same topological properties. In particular, we explore whether cost-integration answers 
this problem. Specifically, we consider whether cost-integration is a useful way of disentangling 
weighted cost from topology. 

The paper is organized as follows. We first introduce, in section 2, some of the notation and 
basic concepts that will be used throughout the paper. In section 3, we describe the two general 
families of topological measures for weighted networks, which are (i) the weighted and (ii) cost- 
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integrated metrics. These quantities are first defined for a single network. The main contribution of 
this paper is then reported in section 4, where different approaches to weighted network comparison 
are outlined, using theoretical results and simple examples. Section 5 describes an application of 
these techniques to a repeated measures fMRI task investigating working memory. This also allows 
us to illustrate a Monte Carlo (MC) sampling scheme to approximate the different measures of 
interest. In section 6, we discuss the findings of this paper in light of the current utilization of 
networks in the biomedical sciences. Finally, we close with a set of recommendations on how 
to conduct weighted network analysis in practice and how to report the findings arising from 
this type of research. An R package entitled NetworkAnalysis (http://CRAN.R-project.org/ 
package=NetworkAnalysis) has been developed that makes available the methods discussed in 
this paper. 

2 Network Types and Topologies 

2.1 Unweighted, Weighted and Fully Weighted Networks 

For clarity of exposition and consistency with the previous literature, we will here employ the no- 
tation used by Kolaczyk (2009). A comprehensive introduction to the theory of complex networks 
can be found in Newman (2010). In the following, the terms metrics and measures will be used 
interchangeably to refer to a function quantifying the topological structure of a network. Our 
use of the terms metric and measure is unrelated to the mathematical definitions of these con- 
cepts in topology and measure theory, respectively. Similarly, we here utilize the graph-theoretical 
definition of the term cost, which is not related to its use in a probabilistic setting. 

An unweighted undirected graph or network G is formally defined as an ordered pair (V,£), 
where V is a set of vertices, points or nodes, and £" is a set of edges or connections linking pairs 
of nodes. Therefore £ C V (S) V, where (E> is the Cartesian product. The cardinality -i.e. the 
number of elements- of V and £ will be referred to as Ny := \V\ and Ne ■— \£\, respectively, 
where | • | denotes the number of elements in a set, and :— that the left-hand-side is defined as the 
right-hand-side. Moreover, the terms network and graph will be used interchangeably. A graph 
with the maximal number of edges is referred to as a complete or saturated graph. For a given 
network G, we denote the corresponding saturated graph as G^***. The cardinality of the edge set 
of G^'^* is denoted by Nj to distinguish it from Ne- Here, the set 1(G), for any graph G is the 
set of indices of all possible edges in G. That is, 

AG) -.^{[i,]) ■.l<i<j <Nv}. (1) 

This notation for the set of indices of all possible edges in G will be useful when describing the 
topology of G based on its shortest paths. 

Weighted undirected graphs will be denoted by the triple G := (V, W), where W{G) is a set 
of weights, whose elements are indexed by the entries in £{G), such that 

Wei = WvjVh, (2) 

for some edge := VjVh- Thus, every weighted undirected graph will necessarily satisfy 

Ne^Nw< Ni, (3) 

where Nw ■= |yV(G)|. The weight set populates a symmetric matrix W, whose diagonal elements 
are null. Graphs that satisfy Nw = Nj will be referred to as fully weighted graphs. Note that, in 
general, we will not draw an explicit difference between a weighted and an unweighted network 
through our notation. However, which one we are referring to should be understandable from the 
context. 

There are a wide range of different weighted measures of internodal association. Our method- 
ological development, in this paper, applies to any choice of association metric. This includes 
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correlation coefficients, partial correlations, synchronization likelihoods and others. For simplic- 
ity, we will assume that the association weights Wj^ 's lie in the unit interval, [0, 1]. Roughly, these 
standardized weights, Wij, can be interpreted as the strength of the association between nodes i 
and j, with larger values indicating a greater level of association. Such standardization can be 
obtained straightforwardly, in practice. For the case of the Pearson's correlation coefficient Vij, 
for example, the standardized weights can be defined as. 



Note that the use of such a standardization of correlation coefficients potentially leads to two major 
pitfalls. Firstly, since negative correlations arc; transformed into positive measures of association, 
it follows that we are amalgamating different subsets of edges, which may play very different 
roles. That is, while subnetworks of negatively correlated vertices may reflect inhibitory processes, 
subnetworks of positively correlated vertices may reflect excitatory processes. Secondly, since pairs 
of vertices linked by a small amount of correlation, either positive or negative, will be transformed 
to take a value close to 0.5; it follows that we may be introducing a spurious amount of random 
noise in such a weighted network analysis, as correlation coefficients close to zero are likely to be 
non-signiflcant. Our approach to weighted network analysis in this paper, however, centres on 
thresholding the weighted networks of interest and therefore does not explicitly take into account 
the direction of the association. Moreover, our focus will be on fully weighted networks, such as a 
standardized correlation matrix, where all entries are greater than 0. Therefore, in the sequel, G 
will refer to a fully weighted graph, except when specified otherwise. We will discuss the use and 
limitations of cost-integration for non-fully weighted graphs in the discussion. 

2.2 Classical Measures of Network Topology 

A wide range of network topological metrics have been proposed in the literature (see Rubinov 
and Sporns, 2010, for a review). Two types of measures are generally of interest, which are 
sometimes referred to as (i) integration metrics and (ii) specialization metrics. The former category 
of topological measures quantifies a network's capacity to transfer information globally, whereas 
the latter reflects a network's capacity to transfer information locally. This distinction originated 
with the work of Watts and Strogatz (1998), who considered the characteristic path length (CPL), 
on one hand, and the clustering coefficient (CC), on the other hand, as measures of global and 
local information transfer, respectively. Although these metrics have been successfully used in 
a wide range of settings. Latora and Marchiori (2001) have introduced two analog metrics: the 
global and local efficiencies, which will be more useful in our context. These two measures retain 
the interpretation of the CPL and CC, while being applicable to a wider range of networks. 
Specifically, global and local efficiency metrics can be computed for any network, irrespective of 
their level of sparsity, which is not true for CPL and CC. (That is, CPL becomes infinite when a 
graph is disconnected and CC becomes undefined when a vertex has no neighbors -that is, when 
a node is isolated.) Throughout this paper, and following other authors (Achard and Bullmorc, 
2007), we will therefore focus on the comparison of families of networks, whose topologies are 
characterized by efficiency metrics. 

One of the remarkable aspects of global and local efficiencies is that they can both be subsumed 
under the general concept of information transfer efficiency, which is defined for any unweighted 
graph G = {V,£) -connected or disconnected- as (Latora and Marchiori, 2001), 



where the summation over the set T(G) is over all the pairs of indices as in equation (1), 

and dij denotes the length of the shortest path between vertices i and j in the adjacency matrix 
of G, with dij := oo when these two nodes are not connected. The summation over j ^ i includes 




(4) 
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all indices between 1 and Ny different from i. The global and local efficiencies of network G are 
then readily derived from equation (5), such that 

£;G'°(G):=£(G), and E^^^^iG) := —J^ E{Gi), (6) 

i—l 

where Gi is the subgraph of G that includes all the neighbors of the i*^ node. That is, V(Gj) := 
{vj e G\vj ^ Vi}, where Vj ~ Vi signifies that nodes i and j are connected. By convention, we have 
Vi ^ V{Gi) (see Latora and Marchiori, 2001, 2003). Note that both global and local efficiencies are 
normalized quantities with values in the unit interval -that is E{G) £ [0, 1]. The global efficiency 
of a graph G can be interpreted as the average 'speed' of information transfer between any pair of 
nodes in G, with a high value of E'^^°{G) indicating a high average 'speed', and therefore efficient 
information transfer. Similarly, the local efficiency of a graph G can be interpreted as the average 
global efficiency of the Ny subgraphs of G, where again a high value for £^^°'^(G) implies efficient 
local information transfer, on average. 

We have used E{G) to denote the efficiency metric of the unweighted graph G. This should 
be distinguished from the graph-theoretical concept of the edge set, which we have denoted £{G). 
Since both quantities are functions of G, we have emphasized this distinction through our notation. 
Note also that we will make use of the expectation operator from probability theory, which will 
be denoted by E[-]. For simplicity, all our development, examples and technical results will be 
based on the general efficiency described in equation (5). However, these methods could readily be 
extended to both global and local efficiencies. In fact, most our discussion applies to all topological 
metrics that can be computed for any level of sparsity. We will discuss the generalization of our 
results to other topological measures in section 6. 



2.3 Cost and Weighted Cost 

In network analysis, it is often of interest to quantify the cost or wiring cost of an unweighted 
graph. In this section, we extend this concept to weighted networks. This generalized version of 
cost will be termed the weighted cost or weighted density. 

The cost or density, K := K{G), of an unweighted network G = {V,£) quantifies the relative 
number of connections in G as a proportion of the number of edges contained in the A''y-matched 
saturated network G^^*. That is, 

where Nj := Nv{Nv — l)/2. The computation of the cost of a network G implicitly refers to the 
adjacency matrix A of that network. Hence, we can reformulate the definition in equation (7) by 
explicitly using A as follows, 

= ^N^Nt^Zt = iv7 ^ "*^' 
Z^i=l l^jjti "-ij 1(G) 

where the afj'^'s denote the elements of the adjacency matrix A^"**, which represents a saturated 
network on Ny nodes. 

Similarly, it will be of interest to quantify the cost of a weighted network, which will be referred 

to as Kw{G). We define it by generalizing the relationship between an unweighted graph and its 
adjacency matrix in order to apply it to weighted graphs and their association matrices. However, 
to extend the concept of cost to a real-valued association matrix, say W, we need to formalize 
what we mean by a saturated weighted graph. A natural choice is to define W^**' as a matrix of 
order (Ny x Ny) with unit entries. Formally, W^*** := I(NvxNv)- Using this saturated association 
matrix, we can now define the cost of a weighted graph as follows, 

^^(^) ^ATv^^^ Sat = TV, 2^ "''^■' 
Z^i=l l^jjii ^ij 1(G) 
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where wfj**'s are the elements of W . The non-standardized version of the cost of a weighted 
network in equation (9) was introduced by De Vico Fallani et al. (2008). Thus, the weighted cost 
of G = (V, 5, W) is the mean of the off-diagonal elements in W, populated by the set W. This is 
reminiscent of our starting point in equation (7), where the same observation can be made about 
unweighted networks. In the sequel, the concept of weighted cost will be used interchangeably with 
the phrase connectivity strength. Note that depending upon which standardization one chooses, 
one may obtain different types of weighted costs. In particular, Kyy could also be standardized 
with respect to the number of elements in W. This would produce a different measure. In 
this paper, we will assume that the networks under consideration are fully weighted, such that 
Ne = Nw = Nj, and therefore these two types of weighted costs are equivalent. 



3 Measures of Weighted Network Topology 

There is currently no guidance in the literature on how to quantify the topological aspects of a 
weighted network. We review here two approaches to this problem: (i) weighted, and (ii) cost- 
integrated metrics of network topology. We describe and define these two families of measures, in 
turn. 



3.1 Weighted Measures 

A natural approach to the problem of quantifying the topology of weighted networks is to translate 
unweighted measures, such as efficiency metrics, for example, into a weighted format. This is a very 
general procedure, which has been introdticcd by several authors including Latora and Marchiori 
(2001) and Rubinov and Sporns (2010). Weighted versions of classical metrics commonly rely on 
the definition of a weighted shortest path. For unweighted networks, the shortest path dij between 
nodes i and j in G = (V, £) is defined as the following minimization, 

where Pij (G) is the set of all paths between nodes i and j that are subgraphs of G. A subgraph 
Pij C G is a path if and only ii i,j G V{Pij) such that 

SiPij) = {ia,ab,...,yz,zj} , (11) 

where each pair of letters stands for an edge. One can similarly define a weighted shortest path, 
dj^ , for some weighted graph G = (V, f , W) as follows, 

where the weighted edge set of a path now takes the form, 

W(Pij) = {Wia, Wab, Wy^^W^i] , (13) 

using the notational convention introduced in equation (2). Since we have normalized the associa- 
tion weights, uijj's, the real-valued function /(•) is restricted to a map of the form / : [0, 1] ^ [0, 1]. 
A convenient choice of /(•) is the inverse function, f{wij) := l/wij. It now suffices to use our 
chosen definition of the weighted shortest path d^ , in order to obtain a weighted version of the 
general efficiency metric in equation (5), which gives 



Note that weighted efficiency is not necessarily bounded between and 1. Here, we have d^ £ 
[min(wjj),oo], regardless of whether or not the Wjj's have been standardized. The case dij := oo 
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may occur when there docs not exist a path between i and j. However, since we have restricted the 
scope of this paper to fully weighted matrices, where Wij € (0, 1) holds for every pair of indices, it 
follows that Ew e 1R+ := (0, oo). 

3.2 Cost-integrated Mecisures 

A second approach to the problem of quantifying the topology of weighted networks proceeds 
by integrating the metric of interest with respect to cost. Here, some authors have integrated 
over a subset of the cost range (see Achard and Bullmore, 2007, for example), whereas others 
have integrated over the entire cost domain (He et al., 2G09a). This second family of metrics will 
be referred to as cost-integrated measures. Given a weighted graph G = {V,£,yV), the general 
efficiency of equation (5) can be defined as follows, 

Ek{G):= J E{^{G,k))p{k)dk, (15) 

nK(G) 

where cost is treated as a discrete random variable K . with realizations in lower case, and p{k) 
denotes the probability density function of K. Since K is discrete, it can only take a countably 
finite number of values, which is the following set, 

where, as before, Ni := Nv{Nv — l)/2 = \Q.k\- It will be useful to treat VLk as an ordered set, k, 
whose elements, fci's, are arranged in increasing order and indexed by i = 1, . . . , iV/. 

The function ^{G, k) in equation (15) is a thresholding function, which takes a weighted undi- 
rected network and a level of wiring cost as arguments, and returns an unweighted network. We 
defer a full discussion of 7 to Appendix A, where we describe its definition in more detail. This 
function is based on the percentile ranks of the elements of W, where tied ranks are resolved by 
assigning the corresponding ordering of the elements' indices. Since there is no prior knowledge 
about which values of K should be favored, we specify a uniform distribution over fix- In equation 
(16), we have excluded the null cost for standardization purposes. Since any edge-based topology 
of interest will be zero when K = Q, this particular value is irrelevant when comparing different 
populations of networks. In example 3, we will also sec that this exclusion of the point mass at 
A' = ensures a more satisfying standardization of Ek(G). 

Since K is treated as a discrete random variable, we can define its probability mass function. 
As no particular cost levels are favored, K is given a discrete uniform distribution, such that 

K ~ DisUnif(Oif ), (17) 

where each element of VLk has an identical probability of occurrence, which, in our case, is equiv- 
alent to 

for every A; G VLk- The; theoretical integration in equation (15) is therefore a weighted summation 
over a finite set of atoms (see Billingsley, 1995), and may be computed as follows, 

Ek{G) = J2 E{l{G, hMh) = ^ E ^*))- (19) 

t=i ^ t=i 

where the index t runs over the elements of CIk described in (16). 

More generally, cost-integrated metrics can be defined with respect to a subset of the cost 
regimen. This perspective on the problem of weighted network comparison has been utilized by 
several authors (Eguiluz et al., 2005, Achard and Bullmore, 2007, Supekar et al., 2009). In our 
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notation, a subset of the cost levels will be indicated by an interval of the form C fi/^, 

which refers to a finite number of values of K, satisfying k- < k < k+. Integration over that 
subset is then defined as 

E[k_,k+]{G)^ j E{j{G,k))p{k\k-,k+)dk, (20) 

where the probability mass function on K is normalized with respect to the chosen domain of 
integration [k-,k+], such that p{k\k-,k+) = l/{Nik+ — Nik- + 1), for every k in that interval. 
The computational formula for this generalization of equation (19) is then given by 

+ t=Nik- 

which follows from Njki = I, using the definition of cost in equation (7). Note that the value of 

the conditional probability p{k\k^, k^) will be different if semi-open intervals such as (fc_, fc+] are 
considered, instead of closed ones. This is due to the fact that the interval of interest is over a set 
of discrete values, as opposed to a subset of the real line. As a special case, this notation can also 
handle the estimation of a particular topological metric at a single cost level, say fco- In such cases, 
the interval of interest becomes [ko, ko]. Our notation makes explicit the fact that integration over 
a subset of the full cost regimen, is conditional on the choice of such a subset. 

Since K has been treated as a random variable and because E{'j{G, K)) is a fimction of 
K, it follows that E{'-f{G, Kj) is also a random variable. The integral Ek{G) can therefore be 
seen as the expectation of E{'y{G,K)) with respect to the distribution of K. This probabilistic 
treatment of cost-integrated metrics will be particularly helpful when considering how to estimate 
these quantities, as a Monte Carlo (MC) sampling scheme can readily be devised in order to 
approximate Ek{G), when the network of interest is too large to be computed exactly. More 
details about this sampling scheme are given in Appendix A. 

4 Pros and Cons of Integrating over Cost Levels 

We now turn to the main question tackled in this paper: Is it useful to integrate over the different 
cost levels of a particular weighted network? In order to answer this question, we briefly consider 
some of the alternatives to this approach. This consists of (i) fixing a cutoff point, (ii) fixing a cost 
regimen, (iii) integrating over all cost levels, and (iv) directly using weighted topological metrics. 
Our comparison of these four approaches is substantiated by some simple examples, synthetic 
data sets, and theoretical results. For convenience, we will solely treat the case of two weighted 
networks in this section. Extensions of these ideas to the case of several populations of networks 
will be discussed in section 6. 

4.1 Fixing a Cutoff Threshold 

The simplest way of comparing the topology of weighted networks is to threshold the corresponding 
association matrices at a specific value, and evaluate the resulting discrete topologies. It is instruc- 
tive to study the consequences of such a naive thresholding on two networks with proportional 
association matrices, as we describe in the following example. 

Example 1. Let two weighted networks Gi = (V, £, Wi) and G2 = (V, £, W2), with standardized 
association matrices denoted Wi and W2, respectively; such that every Wij^k S (0,1) where 
k = 1,2 labels the two graphs under scrutiny. In addition, assume that 

Wi:=aW2, (22) 

where a G (0, 1) is a scalar. That is, the association matrix of G2 is simply proportional to that of 
Gi. Two such association matrices have been discussed in the introduction and were illustrated 
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in panel (a) of Figure 1. Note that the relationship in equation (22) implies that the diagonal 
elements of Wi are not standardized to 1.0. However, the topology and cost of weighted networks 
solely depend on the off-diagonal elements of such association matrices. Therefore, differences 
in the diagonal elements do not pertain to this discussion. Interestingly, it is easy to show that 
proportionality in association matrices implies proportionality in weighted cost. Using equation 
(9), we have 

Kw{G2) = ^ ""'^J.i = ^Kw{Gi), (23) 

^ X(Gi) 

since a is applied elementwise. Therefore, Kyir{G2) > Ky[r{Gi) as by assumption < a < 1. 

A naive approach to the problem of comparing the topologies of these two networks may proceed 
by thresholding Wi and W2 at a particular value, say c*, as was done in the introduction. If 
we compare these networks in terms of global efficiency, straightforward computation of the two 
corresponding quantities shows that we necessarily have 

E{k{G2,c*))>E{k{Guc*)), (24) 

for any c* G [0, 1], where k(G/c, c*) := I{Wfc > c*}. This follows since G2, thresholded at c* has 
all the edges of k{Gi, c*), as well as additional links owing to its weighted cost being higher. The 
relationship in equation (24) is then deduced from the monotonicity of the efficiency function with 
respect to cost. Note that these inequalities would hold for both local and global efficiencies, or any 
other topological metric, which is a monotonic increasing function of the cost level. Therefore, 
example 1 has shown that thresholding weighted graphs at a fixed cutoff point is misleading, 
since graphs with higher weighted cost will tend to be classified as having higher levels of global 
efficiency. This problem can be remedied by fixing cost levels instead of cutoff points. 

4.2 Fixing a Cost Level 

A natural approach to the problem of separating cost from topology is to choose a particular 
cost level. This may be a single value or a subset of the cost regimen. Such a strategy has been 
adopted by several authors (see Eguiluz et al., 2005, Achard and Bullmore, 2007, Supekar et al., 
2009, for examples). One of the original justifications for conditioning over a subset of the cost 
regimen was that topological metrics such as CPL or CC cannot be computed for disconnected 
networks, thereby making it impossible to calculate these quantities for small cost levels. However, 
since comparable global and local topological properties can also be measured using the efficiency 
metrics introduced by Latora and Marchiori (2001), such problems do not arise when using these 
topological metrics. We illustrate the consequences of integrating over a subset of the range of K 
with a real data example, where the original data has been transformed. We have constructed a 
pathological case, which shows that integrating over a subset of the cost levels can fail to distinguish 
between topologically distinct weighted networks. 

Example 2. We here consider a single functional connectivity matrix W, corresponding to the 
mean statistical parametric network (SPN) of a previously published data set, as described in 
section 5 (Ginestet and Simmons, 2011). The matrix W was transformed in order to produce 
two other matrices with either a regular or a hybrid structure, denoted by Wjog := Fi.cg(W) and 
Whyb '■— ■fhyb(W), respectively. The functions i^reg and i^hyb simply re-organize the position 
of the entries in W, as can be seen from Figure 2. The choice of these transformations was 
constrained by the following prescriptions, 

7(aeg,fc') =7(Ghyb,fc') and 7 (G,eg, fc") = 7 (Ghyb, fc") , (25) 

for cost levels k' := 0.25 and k" := 0.75, respectively. That is, the adjacency matrices corre- 
sponding to costs k' and k" are identical for Wrog and Whyb- The effect of the functions F^cg 
and i^hyb was to create different layers of topological structures that vary according to wiring 
cost. The hybrid matrix was composed of alternating layers of random and regular topologies. 
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\ 



Figure 2. Simulation framework for the counterexample in section 4.2. The small-world correlation 
matrix W is transformed into a regular and a hybrid matrix, denoted Wrog and Whyb- The regular 
matrix exhibits a lattice-like topology throughout its cost range, whereas Whyb consists of alternating 
topological layers of random and regular structures. The entries in both matrices have been arranged in 
decreasing order from the diagonal, to facilitate visualization. 



Roughly, the three layers of the hybrid network corresponding to an hybrid association matrix can 
approximately be described as follows, 



{Random if fc G [0, fc'], 
Regular if fc e (fc', fc"]. 
Random if fc G (fc", 1.0]; 



(26) 



for every fc e [0,1], where fc can only take a finite number of values in the unit interval. The 
regular matrix, by contrast, was built as three layers of regular topologies. That is. 



topology (7(G'rog, k)) = Regular, 



(27) 



for every fc G [0, 1]. The random and regular layers were constructed in a standard fashion (see 
Ginestet and Simmons, 2011). Matrices W, Wj-cg and Whyb corresponding to weight sets W, 
Wreg and Whyb, are represented in Figure 2 with the corresponding adjacency matrices resulting 
from different choices of cost levels. 

By construction, the weighted graphs Greg = (V,f,Wreg) and Ghyb = (V,f, Whyb) have iden- 
tical levels of general efficiency for the cost levels comprised in the interval [fc',fc"]. Therefore, 
integrating over that interval gives the same result for both graphs: 



[k',k' 



AG, 



rog/ 



E[k' .k"\{Ghyh) — 0.708, 



(28) 



where = means approximately. By contrast, the general efficiencies of these two networks differ 
substantially when integrating over the full range of cost, i.e. [0, 1]. This gives 



EKiGrcg) = 0.662 and £;i^(Ghyb) ^ 0.679. 



(29) 



This is as expected, since the hybrid network has several layers of random topologies, which renders 
it more globally efficient than Greg- 



Example 2 illustrates the problems associated with integrating over a subset of the cost regimen. 
By doing so, we are potentially omitting substantial topological differences between the networks 
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of interest at other cost levels. The difference in Ek between Greg and Ghyb reported in that 
counterexample may not appear very large. However, these two networks could have represented 
the mean networks of two populations of interest. Providing that the pool of subjects is sufficiently 
large, such topological differences could be found to be statistically significant. By contrast, 
comparison of these two networks on the basis of the full cost regimen yielded answers, which were 
exactly identical, thus nullifying any statistical test of group differences. Naturally, this example 
could have been constructed in the opposite direction in order to show that networks that seem 
to differ topologically for some cost subsets are, in fact, identical when integrating over the full 
cost regimen. 

Fixing a cost level or a subset of the cost regimen therefore suffers from two main problems. 
Firstly, the arbitrariness of the choice of a specific cost subset will generally be difficult to jus- 
tify from either a theoretical or a practical perspective. Secondly, as we have illustrated with 
example 2, considering only a subset of the cost potentially omits topological differences, which 
are solely visible at other cost levels. Thus, any network analysis using this strategy can only 
draw conclusions that are conditional on the choice of cost subset, and this dependence should be 
made explicit when reporting the results of such analyses. Nonetheless, fix;ing a particular cost 
subset successfully satisfies one of our desiderata, which was to disentangle differences in cost from 
differences in topology. That is, weighted networks' topologies can be compared irrespective of 
cost differences, by conditioning on some subset of the cost levels. This invariance property will 
be made mathematically more precise in the next section. 



4.3 Integrating over Cost levels 

From a statistical perspective, the problem of isolating topology from connectivity strength may 

be rcformiilated as evaluating topological differences while 'controlling' for cost, where these two 
quantities are treated as random variables. A natural starting point is to consider weighted 
networks whose association matrices are proportional to each other, as in the ensuing example. 

Example 3. As a simple example, consider the following problem. Let two weighted networks 
Gi = {V,S,Wi) and G2 = (V,£,W2), be characterized by the following standardized association 
matrices: 

0.0 Wi2,2 



Wi 



0.0 Wi2,l 

0.0 



and W2 := 



0.0 



(30) 



where we assume that u>i2,i and ^12,2 are comprised in the open interval (0,1). Here, there 
are only two levels of cost, K e {0, 1}. Trivially, Gi and G2 can therefore be shown to exhibit 
identical general efficiency for these two cost levels. Since our proposed formula for cost-integrated 
topological measures in equation (19) does not include the null cost, we simply have Qk = {1}, 
which implies that both graphs attain the maximal level of cost-integrated efficiency. That is, 

Ek{Gi)=Ek{G2) = 1.0. (31) 

This simple example serves as a justification for our exclusion of the null cost from the set in 
equation (16). Including the null cost would result in = 0.5 for these two basic networks, which 
docs not appear satisfying. Crucially, the equality in (31) does not depend on the relationship 
between Wi2,i and ^12,2- That is, differences in weighted cost have no impact on cost-integrated 
topology. We now return to the case studied in example 1 in order to elucidate the exact effect of 
cost-integration. 

Example 1 (Continued). In this example, we considered two networks with proportional associ- 
ation matrices, satisfying Wi := q;W2. An application of the cost-integrated metrics described in 
equation (19) to the networks of this example gives the following equalities, 

Ek{Wi) = Ek{qW2) = Ek{W2). (32) 

That is, when integrating with respect to the cost levels, we are evaluating the efficiencies of Gi 
and G2 at Nj discrete points. At each of these points, the efficiency of the two networks will be 
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identical, because Wi is proportional to W2 and therefore the same sets of edges will be selected. 
Thus, Gi and G2 have identical cost-integrated efficiencies. 

The equalities derived in these two examples can be shown to hold in a more general sense. 
The invariance of cost-integrated efhciency turns out to be true for any monotonic (increasing or 
decreasing) transformation of the association matrix and applies to any topological metric that 
takes an unweighted graph as an argument, as formally stated in the following result. 

Proposition 1. Let a weighted undirected graph G = (V,f ,yV). For any monotonic function 
h{-) acting elementwise on a real-valued matrix W, corresponding to the weight set W, and any 
topological metric E, the cost-integrated version of that metric, denoted Ek, satisfies 

Ek{W) = EK{h{W)), (33) 

where we have used the association matrix, W, as a proxy notation for graph G. 

A proof of this result is provided in Appendix B. It relies on the idea that the evaluation 
of a weighted network solely depends on the ranking of the off-diagonal elements of W (i.e. the 
ranking of the elements in W), and that the ranks of a set of values are independent of a monotonic 
transformation of these values. Note that the arguments used in Appendix B do not rely on the 
definition of E. Therefore, proposition 1 is true for any cost-integrated topological metrics -i.e. a 
metric originally defined in a discrete setting for an unweighted graph, and integrated with respect 
to cost, when applied to a weighted network. Note also that proposition 1 only holds for any level 
of sparsity in G if the thresholding function 7(6*, k) used in the computation of a cost-integrated 
metric preserves the original ordering of elements in W with tied ranks, using their indices. In 
general, however, sparse networks may better be dealt with, in this context, by adjusting the size 
of the integration domain. 

Proposition 1 encapsulates both the advantages and limitations of cost-integrated topological 
metrics. Two weighted networks, whose topologies are roughly identical at every cost level will be 
given identical scores under this family of metrics, irrespective of cost differences. Cost-integrated 
metrics are therefore successful at winnowing topology from connectivity strength. Another sin- 
gular advantage of this approach is that we obtain a measure, which is invariant under any nor- 
malization or standardization of the original data. That is, any functions that simply rescale or 
shift the association weights, in order to ensure that they are comprised in the unit interval, for 
instance, will have no eff'ect on the value of the cost-integrated topological measures. 

However, proposition 1 also demonstrates the limitation of such an approach. One can easily 
see that such cost-integration will potentially mask some cost-specific topological differences, as 
illustrated in example 2. In addition, when cost-integrated topological metrics are used for network 
comparison, this requires that the sizes of the weight sets of different networks are identical. 
Similarly, the presence of multiplicities in the ranks of the weights may also cause problems, as 
this would artificially induce a random topological structure, since weights with equal ranks would 
be randomly allocated to different cost levels. We will further discuss these limitations in the 
conclusion of this paper. 



4.4 Using a Weighted Metric 

A seemingly natural way of amalgamating connectivity strength and topological characteristics is 

by directly considering weighted topological metrics, such as the weighted global efficiency, Ew, 
introduced in equation (14). Unfortunately, we here prove that such an approach suffers from 
a serious limitation, which could potentially dissuade researchers from using this particular type 
of metrics. With the next proposition, we show that in a wide range of settings, the weighted 
efficiency is simply equivalent to the weighted cost of the graph of interest. 

Proposition 2. For any weighted graph G = {V,S,W), whose weight set is denoted by W{G), if 
we have ^ 

min Wij > - max Wij, (34) 
ju,jeW(G) 2 tu,jeW(G) 
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Sagittal SPNj 

Sup. Sup. Sup. Sup. 




Inf. Inf. Inf. Inf. 



0-back 1-back 2-back 3-back 

Figure 3. Mean Statistical Parametric Networks (SPNj) over the 4 levels of the A*'-back task, in the 
sagittal plane, based on wavelet coefficients in the 0.01-0.03Hz frequency band, with FDR correction (base 
rate ao ~ .05). Locations of the nodes correspond to the stereotaxic centroids of the corresponding cortical 
regions. The inferior-superior orientation axis is indicated in italics. The size of each node is proportional 
to its degree. 

then 

Ew{G)=Kw{G). (35) 

This result can be proved by contradiction, as demonstrated in Appendix E. The hypothesis 
in proposition 2 may at first appear relatively stringent. However, it will encompass a wide 
range of experimental situations. For the real data set described in example 2, the difference 
between maxwij and 2minwij is close to, but not exactly zero. However, we nonetheless have 
Ew — Kw, for that example. Thus, the added value of using the weighted efficiency will, in 
general, be questionable since there exists a strong relationship between this topological measure 
and a simple average of the edge weights. 

These theoretical results and associated counterexamples have therefore highlighted the limi- 
tations of various approaches to the problem of disentangling differences in cost from differences 
in topology. As a result, when comparing several populations of networks, we recommend the 
reporting of differences in weighted costs and differences in cost-integrated topological measures. 
We illustrate this approach with a re-analysis of a previously published fMRI data set. 

5 A/^-back Working Memory Data Set 

In this section, we illustrate our theoretical results with a previously analyzed data set of a 
working memory task based on functional Magnetic Resonance Imaging (fMRI) data (Ginestet 
and Simmons, 2011). In particular, we use this data set for testing our proposed MC sampling 
procedure and for comparing a graph's weighted cost with its cost-integrated and weighted global 
efficiencies. 

5.1 Description 

Ginestet and Simmons (2011) considered topological changes in functional brain networks under 
different levels of cognitive load. Here, we solely give a cursory description of the experimental 
procedure used in this study and refer the reader to the original paper for the full technical 
details. Ginestet and Simmons (2011) constructed networks on the basis of fMRI data gathered 
from 43 healthy adults undergoing a working memory task known as the A'^-back paradigm. In 
this experiment, subjects were shown one letter every two seconds, and were asked to monitor the 
stimuli, in order to indicate by the push of a button whether the current letter was identical to 
the one presented A^ trials previously, where A^ = {1, 2, 3}. A control or null condition was also 
included, the 0-back task, which consisted of simply indicating whether the current letter was an 
X. In this experiment, the subject-specific fMRI images were parcellated into 90 regions of interest 
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Figure 4. Running means of Monte Carlo (MC) estimates for cost-integrated global and local efficiencies 
in panel (a) and (b), respectively, for the 3-back network described in example 2. The grey ribbon 
represents the variability of these estimators at each m = 1, . . . , 5000, using twice the MC standard error. 



That is, ± 2a'^^' for both global and local efficiencies. The dashed lines indicate the exact value of 

Ek- See appendix A for details. 



Am) 



using the Anatomical Automatic Labelling (AAL) template (Tzourio-Mazoyer et al., 2002). The 
BOLD time series were averaged for each AAL region. These regional mean time series were then 
wavelet decomposed. Wavelet coefficients in the low frequency range (0.01-0.03Hz) were selected 
for the main network analysis (see also Achard et al., 2006, for a similar analysis of fMRI data). 
Since the A^-back paradigm contains four experimental levels, we decomposed these time series 
into blocks corresponding to each A^-back condition. As each condition was repeated more than 
once, these blocks were then concatenated. Note that this sequence of processing steps involving 
wavelet decomposition immediately followed by block concatenation was studied by Ginestet and 
Simmons (2011) using simulated data, and was not found to bias the results of the final network 
analysis. 

Vertices in these subject-specific functional networks were chosen to be the 90 AAL regions, 
and the edges were constructed by computing pairwise correlations between each condition-specific 
time series of wavelet coefficients. The results of this construction can be summarized using 
Statistical Parametric Networks (SPNs), as illustrated in Figure 3 (see Ginestet and Simmons, 
2011, for details). SPNs are estimated using a mass- univariate approach, where the edges in a 
population of subject-specific networks are tested for significance using a mixed-effects model, 
and then thresholded using the false discovery rate (Benjamini and Hochberg, 1995, Nichols and 
Hayasaka, 2003). SPNs can be constructed using functions made freely available through the 
R package NetworkAnalysis (http : //CRAN .R-project . org/package=NetworkAnalysis). From 
Figure 3, one can observe that the connectivity strength (i.e. weighted cost or averaged correlation 
coefficient) of the functional networks in each condition tend to diminish as subjects experience 
greater cognitive load. 

5.2 Monte Carlo (MC) Estimation 

A full description of the theory supporting MC estimation in this context is provided in Appendix 
A. MC techniques are here used to speed up the computation required when estimating our 
proposed cost-integrated measures. Figure 4 shows the convergence of to Ex, for a medium- 
sized weighted network derived from fMRI data on the working memory task described in example 
2. The results are provided for both global and local efficiencies. Each plot in Figure 4 shows the 
running mean plus or minus twice the running MC standard error, which are defined for the cost- 
integrated efficiencies, as i?^'"* and (w^'')^^^, respectively, where m = 1, . . . , 5000. (See Appendix 
A for details.) In Figure 4, we also report the exact values of Ek using formula (19) by dashed 
lines. 

In all the cases studied, the MC estimates compared favorably with the exact integrals after 
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approximately a quarter of the number of computations rc;quircd for the exact calculations. That 
is, the exact derivation of necessitates Nj = 4005 evaluations of the global or local efficiency. 
By contrast, MC estimates based on about 1000 samples appear to provide reasonably good 
approximations of these quantities, as indicated by the small MC standard error. This constitutes 
a non-negligible computational gain. The MC standard error, which is derived as a by-product of 
these computations could then be used as an indicator of the uncertainty associated with these 
estimates in a Bayesian hierarchical model, where uncertainty is propagated from the data to the 
population's parameters of interest. 

A simple alternative to MC averaging, in our context, would be to construct a mesh of the unit 
interval and to approximate the desired integral by a weighted sum of the values of the topological 
metric of interest at the midpoints of that mesh. The latter method is generally referred to as 
the Gauss-Kronrod quadrature formula (see chapter 2 of Minka, 2001, for a review of integration 
methods). While this method is very efficient for simple functions, it becomes rapidly unwieldy for 
complex ones, as it requires an increasingly refined mesh to ensure good interpolation. Moreover, 
since the Gauss-Kronrod is a deterministic algorithm, it does not provide a measure of the accuracy 
of the estimation. By contrast, a MC approach ensures asymptotic convergence for any level of 
complexity and also produces precise confidence bands. (See Appendix A for details.) 



5.3 Evaluation and Comparison 



Following the statistical framework used in the original analysis of this data set (Ginestet and 

Simmons, 2011), we tested for the statistical significance of the TV-back factor on different topo- 
logical metrics using a mixed-effects model. We here have n = 43 subjects and J = 4 experimental 
conditions. Using the formalism introduced by Laird and Ware (1982), we have 



y, = X,/3 + Z,b, . 



l,...,n, 
7V(0,a2l), 



(36) 



where 



[yn 



, yij]"^ is a subject-specific vector of topological metrics of interest, j3 := 



. . . is a vector of fixed effect, which do not vary over subjects, b, := bn is a subject- 

specific random effect and e, := [en, . . . , Cjj]^ are the residuals. Finally, the matrices Xj's and 
Zj's are given the following specification, 



10 
110 
10 10 
10 1 



and 



(37) 



for every i = 1, . . . ,n. The effect of the TV-back factor was then evaluated using Wald's F-test. All 
these analyses were conducted within the R environment using the lmc4 package (sec www.cran.r- 
project.org and Pinhciro and Bates, 2000). Note that the model used here is slightly simpler than 
the one used in Ginestet and Simmons (2011), as the present mixed-effect model was found to be 
better identified than the growth curve model utilized in the original analysis. 

In Figure 5, we report the cost-integrated global efficiencies for this experiment. For illustrative 
purposes, we have computed these quantities for four different choices of domains of integration. 
The E^^'^^{Gij) were here estimated using 1,000 MC samples for each subject in each A^-back 
condition. In panel (a), one can observe a clear increase of the cost-integrated global efficiencies 
as we increase the size of the domain of integration, due to the monotonicity of global efficiency 
with respect to cost. This is a standard property of global efficiency: as graphs become denser, 
their diameter tends to diminish (Bollobas and Riordan, 2003). In Figure 5, one can also note the 
dependence of the inter-subject variability of the cost-integrated metrics on the chosen domain of 
integration. 

We therefore tested for the effect of the A^-back factor on the topological metrics of interest, 
given different domains of integration, in order to evaluate whether such a choice of domain has 
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Table 1. Statistical inference for the mixed-effects 
model described in equation (36) testing for the effect of 
the A''-back factor on different topological variables. For 
cost-integrated global efficiencies, we have separately 
tested four different domains of integration. 



Outcome Variable Domain 



F-statistic^ p-value 



Kw{G) 
Weighted Cost 



3.59 



0.01 



Cost-integrated 
Cost-integrated 

Cost-integrated 
Cost-integrated 



[0, .25] 
[0, .50] 

[0, .75] 
[0, 1.0] 



0.34 
0.24 

0.40 
1.09 



0.79 
0.86 

0.75 
0.35 



^ Wald F-statistic based on model described in equation (36) . 



a systematic impact on the effect of the experimental factor. These tests are based on the mixed- 
effects model described in equation (36) . and wc have reported the results of these statistical tests 
in Table 1. These results do not indicate that the choice of different domains of integration yield 
a systematic bias in statistical inference. As was reported by Ginestet and Simmons (2011), the 
weighted cost was found to be systematically affected by the A^-back factor (Wald F = 3.59, dfi = 
3,df2 = 126, p = 0.01). However, none of the cost-integrated global efficiencies appeared to 
be significantly influenced by the experimental factor. Most importantly, the use of different 
domains of integration did not seem to affect the results. Integration over the entire cost domain, 
however, resulted in a larger F-statistic, which may be explained by the lower amount of variability 
characterizing cost-integration over larger domains, as can be observed in Figure (5). 

In addition, in Table 1, we have also reported the F-statistic for the effect of the A-back 
factor on the weighted cost. The subject-specific network's weighted costs were found to be 
significantly influenced by the level of the experimental factor, as is immediately visible from the 
mean SPNs reported in Figure 3. The separation of the differences in cost from the differences 
in topology that results from the use of a cost-integrated topological metric is best illustrated by 
the interaction plots in Figure 6, where ensembles of global efficiencies corresponding to different 
costs are represented for the four levels of the experimental factors. Note that, here, wc are 
reporting the efficiency metrics for a single level of cost, not integrated over a subset of the cost 
regimen as was done in Figure 5. This is a visual depiction of the A-back factor that corroborates 
the conclusions reached using cost-integrated topological metrics, which stated that topology, as 
measured by global efficiency, does not significantly vary with the experimental factor. 

6 Discussion 

This paper has investigated the effect of thresholding matrices of correlation coefficients or other 
measures of association for the purpose of producing simple unweighted graphs. On the basis of 
the analysis, examples and counterexamples studied in this paper, we make the following method- 
ological recommendations to researchers intending to compare the topological properties of two or 
more populations of weighted networks. 
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Figure 5. Box plots of cost-integrated global efficiencies of flVIRI iV-back networks for four different 

domains of integration. Tiiesc integrals were estimated using MC approximation over 1,000 samples 
for each of the 43 subjects in each of the four experimental condition. Note that different domains of 
cost-integration do not induce any differences in the effect of the experimental factor. For all choices of 
integration domain, there is no apparent significant differences. 



6.1 Summary and Recommendations 

We here summarize the main findings of this paper: (i) fixing a cutoff threshold is not satisfactory, 
because this is fully determined by differences in connectivity strength, as we have shown in section 
4.1; (ii.) fixing a subset of cost levels is not satisfactory, because this potentially omits topological 
differences at other cost levels (see section 4.2); (iii) integrating over the entire cost regimen 
successfully disentangles connectivity strength from topology up to monotonic transformations. 
Specifically, such metrics are invariant to monotonic transformations of the association weights, 
as described in section 4.3; and (iv.) the weighted topological metrics, such as E\y^ appear to be 
too closely related to weighted costs (see section 4.4). 

From a methodological perspective, wc therefore recommend the following. As a preliminary 
step, it is good practice to standardize the association weights, in order to obtain Wij S [0, 1] for 
all Wij's, with large values of the weights corresponding to strong associations. This may facilitate 
comparison across separate network analyses, and ease the interpretation of the results. Secondly, 
the weighted cost, i.e. connectivity strength, of the networks of interest can then be computed for 
all networks. This is central to the rest of the analysis, and should be conducted systematically. 
Moreover, quantitative differences in connectivity strength per se are informative about the brain 
processes at hand, and their experimental relevance should not be neglected. Thirdly, population 
differences in cost-integrated topological metrics may then be evaluated. This will indicate whether 
the topologies of the populations under scrutiny vary significantly irrespective of their differences 
in connectivity strength. This aspect of network analysis could be regarded as qualitative, as this 
reflects the networks' architectural properties. We now expand and discuss some of the remarks 
that were made in section 4.3. 

6.2 Limitations of Cost-integration 

As with any form of averaging, cost-integration ignores cost-specific topological differences. Net- 
works Gi and G2, in example 1, differ in connectivity strength and these differences may also be 
expressed through their cost-dependent respective topologies. That is, as illustrated in example 2, 
certain graphs may not exhibit the same topological structure at different cost levels, and therefore 
integrating over cost may potentially mask these subtle topological differences. Another potential 
problem with cost-integrated quantities is that they may be expensive computationally. The num- 
ber of possible cost levels increases at rate 0{Ny) with respect to the number of vertices in the 
networks of interest. In Appendix A and in section 5, however, we show how such integrals can 
be estimated through MC sampling, which can substantially diminish the required computations. 
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Another potential pitfall which is not directly visible from proposition 1 is that the use of 
cost-integration for the comparison of several populations of networks requires these networks to 
have the same number of positive weights. That is, to be comparable two networks do not simply 
need to possess the same number of vertices, i.e. |V(G'i)| = |V(G2)|) but also should have the 
have the same number of weights, i.e. |W(Gi)| = |>V(G2)|. In this paper, we have re-analyzed an 
fMRI data set, based on correlation matrices, which produce fully weighted networks, for which 
Nj = for every subjects. However, when such a condition does not hold, we recommend the 
selection of a domain of integration that corresponds to the smallest common denominator. That 
is, := minj=i^...^„ |>V(Gi)|, for a given population of n weighted networks denoted Gj. Thus, 
when considering sparser networks, such as structural brain networks, one may still be able to 
control for differences in cost, by integrating over a subset of the cost regimen, which reflects the 
sparsity of the networks under comparison. 

A similar problem may arise if one or several networks in the population of interest have mul- 
tiplicities, i.e. weights that take identical values. Since cost-integration relies on the ranking of 
weights, it follows that one may need to adjust for such multiplicities, otherwise this can lead to 
spurious generation of random topologies. That is. when the tied ranks are resolved by random 
ordering, the allocation of weights with identical values to specific cost levels is random, and there- 
fore artificially create a random topology for these cost levels. For sparse networks, multiplicities 
arc likely to arise around zero. However, for large non-sparse networks, the occurrence of multiplic- 
ities should be evaluated by counting the number of tied ranks in the distribution of the weights. 
In particular, if the two populations of networks that one wishes to compare differ significantly 
in number of tied ranks, then comparison based on cost-integration will be contaminated by an 
artificial level of random topology. 

Another possible limitation of cost-integration is that by integrating over several cost levels, 
we omit to take into account the dependence between the topologies of the different thresholdcd 
graphs. The topological structure of the unweighted networks created by thresholding the original 
weighted graph share the same edges. Arguably, the cumulative nature of this procedure results 
in emphasizing the importance of the set of edges with the largest weights. Once these edges have 
been included into a thresholded graph, they will be retained for the remaining cost levels. This 
is especially true for the topological metrics that we have studied in this paper, since global and 
local efficiencies arc both monotonic functions of cost. 

Finally, one may also be interested in addressing how differences in cost and differences in 
topology interact. By controlling for monotonic differences in wiring cost, we potentially ignore 
how cost differences may contribute to the topological structure. In the supplementary methods, 
we report a different type of integrated topological metrics, which attempts to combine cost and 
topological differences. In this case, cost-integration was weighted with respect to the distance 
between each pair of weights used for thresholding the graph (sec the supplementary methods 
document). Unfortunately, we have shown that this choice of integration exhibits some undesir- 
able properties, in the sense that it tends to give more importance to very low- and very high-cost 
topologies. Further research will therefore be required to produce more relevant topological func- 
tions of weighted graphs, which provide a better understanding of the interaction between cost 
and topology. 

6.3 Future Work 

Most of this paper has focused on the global efficiency metric. Thus, our conclusions and the 
examples studied will not necessarily apply to other topological measures. However, our main 
result (proposition 1) was proved in a very general setting, which is independent of the particular 
formula of the topological metric of interest. Our general conclusion about the usefulness of 
cost-integration when one wishes to disentangle differences in cost from differences in topology 
is therefore valid for any topological metric defined for an unweighted graph. In addition, we 
note that since most weighted metrics are constructed on the basis of the weighted shortest path 
matrix, one may surmise that our second main theoretical result (proposition 2), may hold in a 
more general setting. However, a proof that the equivalence relationship between the weighted 
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Figure 6. Interaction plots of cost-dependent global efficiencies of fMRI networks with respect to the 
levels of the A^-back factor. We here consider five different costs K G {0.05,0.15,0.35,0.55,0.75}. The 
dashed lines represents the cost-specific global efficiencies of each subject, whereas the plain line represents 
cost-specific global efficiencies averaged over the 43 subjects. The flatness of the lines at each cost levels 
suggests that the experimental factor has little effect on the topological structure of these networks. 



version of a topological metric and the weighted cost, for instance, hold for topological measures 
other than the global efficiency would require further work. 

Thus far, we have only considered flat distributions on the space of network costs. However, 
future methodological developments will be needed in order to consider more sophisticated ap- 
proaches to this problem. In particular, the specification of a prior distribution on K should 
take into account the eff'ect size associated with different values of this random variable. When 
considering correlation coefficients, for instance, it can easily be shown that higher values indicate 
larger effects, and it may therefore be preferable to emphasize network comparisons built upon 
the largest correlation coefficients. This may be implemented by integrating network topological 
metrics with respect to a skewed distribution on K, which puts more weight on sparse networks, 
whose edges are better identified. 

One should note that the use of cost-integration when comparing weighted networks is not akin 
to taking into consideration the multilevel or hierarchical nature of a weighted network. Such a 
structural interpretation of the successive thresholding necessary for such an integration is not nec- 
essary to justify the usefulness of the method in controlling for monotonic differences in weighted 
cost. Since the networks of interest 'exist' as weighted networks, their thresholding remains arti- 
ficial and it is not clear whether one can ascribe any substantive meaning to the resulting family 
of thresholded graphs. Further work will therefore be needed in order to better characterize the 
architecture of the ensemble of thresholded discrete networks subtending a weighted graph. 
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Appendices 

A: Monte Carlo (MC) Sampling 



The cost-integrated quantities introduced in section 3.2 may first appear unwieldy to compute, especially 

when considering large graphs. However, the structure of these integrals allows the construction of a 
straightforward MC sampling scheme. This classical approximation method has the advantage of providing 
both an estimate of the quantity of interest and an estimate of the variance of that estimation. For an 
introductory text to MC techniques, see Gilks et al. (1996), and for a more advanced treatment, Robert 
and Casella (2004). 

In order to apply MC sampling theory, we first observe that our integration problem -that is, the 
computation of Ek- can be re-formulated as an expectation. For convenience, we drop any reference to 
the function 7(G, K), and therefore denote the efficiency metric E{'y{G, K)) as E{K). The cost-integrated 
metric Ek can then be expressed as an expectation of E{K) since, straightforwardly, we have 



Ek 



f E{k)p{k)dk = 'Ep[E{K)], (38) 



where fix is the space of all possible costs for G, with \Qk\ = ( 2^) ~ Nv{Nv — l)/2. The expectation 
in (38) is taken with respect to p, the probability density function of K, and explicit reference to G has 
been omitted. It is natural to consider the use of a sample {fci, . . . , km} from p in order to approximate 
Ek by the following empirical average. 



^ m 

-' = ™E^(^0- (39) 



E' 

The approximation JS^' converges to Ek almost surely, by the Strong Law of Large Numbers. In addition, 
providing that E{K) is square-integrable, the speed of convergence can be evaluated by considering 

the theoretical variance of that estimate, 

Var (e^k^^ " rri / (^^''^ ~ i^i^)])^ pC')^'^- (4°) 

10,1] 

which can be approximated by the following MC variance, 

(-irO =iE(^w-i^r)^ (41) 
1=1 

This quantity is of special interest in MC sampling, as it permits the evaluation of the rate of convergence 
of the estimation. It is generally referred to as the MC standard error. Using Slutsky's theorem, it can 
also be shown that as m — ^ 00, the random variable. 



(42) 



tjK — rLK 



"K 

has the probability density function of a Normal variate centered at zero, with unit variance. MC sampling 
is especially useful when the stochastic function that we wish to integrate -here, denoted E{K)- is complex, 
whereas the random variable with respect to which we integrate can easily be sampled. Most topological 
metrics will be of a complex nature i.e. non-linear. By contrast, both K and C will be straightforward 
to sample, since we have specified uniform distributions for both of them. The theory underlying MC 
sampling is general and can therefore be applied to any type of topological metrics. Care, however, should 
be taken when evaluating the properties of very large networks, where the topology may vary substantially 
from one level of cost to another. When confronted with such large networks, however, the MC standard 
error is still a good indicator of the accuracy of such approximations. 

B: Proof of Proposition 1 

In order to prove proposition (1), we first need to give a formal definition of 7(G, fe), for some given 
weighted network G = (V, f,W). This function relies on the concept of rank, which can be formally 
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defined in our context, as follows 

Ny Ny 

Rij (W) ■■=^J2Y1 ^^""'i <'Wuv}, (43) 

where Rij = 1 implies that Wij is the largest weight in W. Here, we have assumed that there are no ties 
in the ranks of W. When ties occur in practice, we recommend to resolve tied ranks by assigning the 
corresponding ordering of the elements' indices. By contrast, resolving tied ranks using random allocation 

can result in introducing a spurious amount of random topology in the networks of interest. The presence 
of tied ranks, however, will generally be indicative of a high level of sparsity, which is better dealt with 
by restricting the domain of integration. 

Computationally, this definition can be simplified if one only considers the upper off-diagonal elements 
of W and omits the division by 2. For our purpose, this definition will be more convenient. These ranks 
can be standardized in order to derive the percentile ranks, 

^'^■(W):-^, (44) 

where Ni is the number of edges in the saturated version of G. Note that the resulting matrices R and P 
of ranks and percentile ranks, respectively, are both symmetric. A good introduction to order statistics, 

ranks and percentile ranks is provided by Lin et al. (2006). 

The function 7(G, k) can now be given a formal definition using the Py 's, such that 

7(G,fc) :=: 7(W,fc) :=X{P(W) < fc}, (45) 

where the indicator function is applied elementwise to matrix P(W), where W is the similarity matrix of 
G. It can hence be seen that the function 7 prescribes an adjacency matrix A'*'^ with the desired cost. 
This can be verified by computing the cost of the corresponding unweighted network G^*-* = (V ,£^''^)., 
where f is the edge set that populates A'''', obtained after application of the 7 function at k. Provided 
that k € fi/f, as defined in equation (16), we have 

KiG^"' ) = ^ E 4' = ^ E (A*'^^ ) < fc} = fc, (46) 

' 1(G) ' 1(G) 

which can be verified by noting that equation (46) is simply the discrete version of the integration of an 
indicator function of the form, T{x < k}dx = dx = k. Using this notation, the proof of proposition 
1 is now straightforward. This demonstration uses the fact that a monotonic function does not modify 
the ranks of its arguments. 

Proof. Recall that the cost-integrated version of E{G), in its computational form, is given by 

Ni-l 

Ek{W) = — E{j{W,h)). (47) 
' t=i 

To demonstrate that Ek(W) = £if (/i(W)), it therefore suffices to show that 

Sk(7(W, kt)) = £k(7(/i(W), kt)), (48) 

for every kt, which further simplifies to the sole requirement that 7(W, fct) = 7(/i(W),fct), for all t = 
l,...,Ni. Prom the definition of the 7 function introduced in equation (45), we have the following 
relationship, 

f{h{W),kt) = I{P{h{W)) < fc} = X I ^'^'^J^^^^^ < fc} . (49) 
However, one can observe that, since h is applied elementwise, we have 

Ny Ny 

Rij{h{W)) = - J2J2^{H^ii) < HvJuv)} = Rij{W), (50) 

u — l VJ^U 

for any monotonic function h. Note that this argument makes no use of the definition of E. This completes 
the proof. □ 
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C: Proof of Proposition 2 

Proof. Wc prove the result by contradiction. Assume that the conclusion does not hold. That is, Ew 7^ 
Kw- By applying the definitions of Ew and Kw in equations (14) and (9), respectively, we have 

Ew{G) ..= ±J2^^^J2^^^ =■■ K^^G). (51) 

1(G) »J 1(G) 

It therefore suffices to show that dYj 7^ w"-^ for at least one of the weights. The weighted shortest path 
is defined in equation (12) as 



dlT := min w^v ■ (52) 



It follows that ^ w^j^ if and only if there exists a path Pij in l^ij{G), which satisfies 

^ wZl < w'/. (53) 

That is, the path P^j is shorter than the direct connection Wij between the i**^ and vertices. Inequality 
(53) can be sandwiched in the following fashion, 

\£iP'j)\( m&x Wij) < V w-J min Wij) , (54) 

wuvee(Pii) 

where |f (Pjj)] denotes the cardinality of Pij. Inverting the entire inequality then gives 



max Wii > 



V w^J >Wij> mm Wij. (55) 
^ , I wijes(G) 
y^u.eeip:^) J 



\£{Pij)\ v.,,es(G) 
However, we clearly have 

^ max Wij > max Wij > min Wij, (56) 

2 wijes(G) \£iPij)\ wijesiG) t««ef(G) 

which contradicts our hypothesis, and proves the claim. □ 
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